IrrigationUpdate Subroutine

SUBROUTINE IrrigationUpdate &
!
(time, sat, Qriver, dtOut, timeOut)

IMPLICIT NONE

!Arguments with intent(in):
TYPE (DateTime),  INTENT (IN) :: time !!current time
TYPE (grid_real), INTENT (IN) :: sat !!soil saturation (0-1)
TYPE (grid_real), INTENT (IN) :: Qriver  !!river discharge at current time
INTEGER (KIND = short), INTENT(IN) :: dtOut !!time step for output

!Argument with intent (inout):
TYPE (DateTime), INTENT(INOUT) :: timeOut

!Local declarations:
INTEGER (KIND = short) :: i,j,k,r,c
REAL (KIND = float) :: flux !!irrigation flux (m/s)
REAL (KIND = float) :: fluxQ !!irrigation disharge (m3/s)
INTEGER (KIND = short) :: doy
REAL (KIND = float)    :: meanSat
REAL (KIND = float) :: Qdownstream

!----------------------------end of declarations-------------------------------

timeString = time

irrigationFlux = 0.

Qirrigation = 0.

cpQriver = Qriver

DO k = 1, fields % count 
    !check soil saturation and irrigation season
    doy = DayOfYear (time)
    fluxQ = 0.
    flux = 0.
    !meanSat = GetMean (sat, maskInteger = fields % elem (k) % mask )
    IF (doy > fields % elem (k) % doy_start .AND. &
        doy < fields % elem (k) % doy_stop  ) THEN ! .AND. &
        !meanSat < fields % elem (k) % sat_max) THEN
        !compute irrigation flux [m/s]
        r = fields % elem (k) % r
        c = fields % elem (k) % c
        
 
        IF (cpQriver % mat (r,c) > fields % elem (k) % e_flow (doy) ) THEN
           fluxQ = cpQriver % mat (r,c) - fields % elem (k) % e_flow (doy)
           !limit to concessed discharge
           IF (fluxQ > fields % elem (k) % max_discharge (doy) ) THEN 
             fluxQ = fields % elem (k) % max_discharge (doy)
           END IF
           
           !store irrigation discharge
           Qirrigation % mat (r,c) =  Qirrigation % mat (r,c) + fluxQ
           
           !conversion m3/s -> m/s
           flux = fluxQ / fields % elem (k) % area  
           
        ELSE
           fluxQ = 0.
           flux = 0. 
        END IF
        
        fields % elem (k) % fluxQ = fluxQ
        
        cpQriver % mat (r,c) = cpQriver % mat (r,c) - fluxQ
        
        !populate irrigation map 
        DO i = 1, irrigationFlux % idim
          DO j = 1, irrigationFlux % jdim
            IF (fields % elem (k) % mask % mat(i,j) /= &
                fields % elem (k) % mask % nodata) THEN
                  irrigationFlux % mat (i,j) = irrigationFlux % mat (i,j) + &
                                              flux * fields % elem (k) % eta
            END IF
          END DO
        END DO
        
    ELSE
        fields % elem (k) % fluxQ = 0. 
    END IF !doy

END DO

!write diverted and downstream released discharge
IF (time == timeOut) THEN
    WRITE (unitIrrigationDiverted,'(a)', advance = 'no') timeString
    WRITE (unitIrrigationDownstream,'(a)', advance = 'no') timeString
    WRITE (unitIrrigationUpstream,'(a)', advance = 'no') timeString
    
    DO k = 1, fields % count 
   
        r = fields % elem (k) % r
        c = fields % elem (k) % c
        Qdownstream = Qriver % mat (r,c) -  Qirrigation % mat (r,c)
        
        IF ( k < fields % count ) THEN
            WRITE (unitIrrigationDiverted,fmt='(" ",e14.7)', advance = 'no') &
                    fields % elem (k) % fluxQ
            WRITE (unitIrrigationDownstream,fmt='(" ",e14.7)', advance = 'no')&
                    Qdownstream
            WRITE (unitIrrigationUpstream,fmt='(" ",e14.7)', advance = 'no')&
                    Qriver % mat (r,c)
        ELSE
            WRITE (unitIrrigationDiverted,fmt='(" ",e14.7)') &
                    fields % elem (k) % fluxQ
            WRITE (unitIrrigationDownstream,fmt='(" ",e14.7)') Qdownstream
            WRITE (unitIrrigationUpstream,fmt='(" ",e14.7)')  Qriver % mat (r,c)
           
        END IF
    
    END DO
timeOut = timeOut + dtOut

END IF

RETURN

END SUBROUTINE IrrigationUpdate